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Abstract 

This paper presents a detailed study of four far- 
field boundary conditions used in solving the single 
airfoil gust response problem The boundary con- 
ditions examined are the partial Sommerfcld radia- 
tion condition with only radial derivatives, the full 
Sominerfeld radiat ion condition with both r adial and 
tangential derivatives, the Bavliss-Turkel condition 
of order one, and the Hagstrom-Hariharan condition 
of order one. The main objectives of the study were 
to determine which far-field boundary condition was 
most accurate, which condition was least sensitive to 
changes in grid, and which condition was best overall 
in terms of both accuracy and efficiency. Through a 
systematic study of the flat plate gust response prob- 
lem, it was determined that the Hagstrom-Hariharan 
condition was most accurate, t lie Bavliss-Turkel con- 
dition was least sensitive to changes in grid, and 
Bavliss-Turkel was best in te rns of both accuracy 
and efficiency. 

I. Introd in tion 

A fundamental problem in unsteady aerody- 
namics and aeroaeoustics is he single airfoil gust 
response problem. This problem has a number 
of challenging features whi< h make it an ideal 
model problem for CAA code validation. Among 
these features are mean flow distortion of the con- 
verted gust, vortex shedding in the wake, un- 
steady pressure on the airfoil surface, and acous- 
tic waves which radiate to die far field. In the 
recent Third Computational Aeroaeoustics Work- 
shop on Benchmark Problems 1 , a series of bench- 
mark solutions was presented for the gust problem 
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using 4 he GLST3D code 2-4 . Numerical solutions 
were also presented at the Workshop by Hixon 1 , et 
al and Kopriva 1 , et al. Additional solutions have 
been presented by Rasetarinera 5 , et al and WangA 
et al. While good agreement has been demonstrated 
for many of the comparisons, some cases have shown 
substant ial discrepancies. 

In an effort to improve the accuracy of the 
GUST3D solutions, a detailed study of four far- 
field boundary conditions was carried out. The 
boundary conditions examined were t he partial Som- 
inerfeld radiation condition with only radial deriva- 
tives, the full Sominerfeld radiation condition 11 with 
both radial and tangential derivatives, the Bavliss- 
Turkel 7 condition of order one. and the Hagstrom- 
Hariharan 8 condition of order one. The main ob- 
jectives of the study were to determine (i) which 
far-field boundary condit ion (FFBC) was most accu- 
rate, (ii) which FFBC was least sensitive to changes 
in grid, and (iii) which FFBC was best overall in 
terms of bot h accuracy and efficiency. To make these 
assessments, we const ructed a matrix of test cases 
for the flat plate gust problem, and then repeatedly 
tested each FFBC on a large variety of grids for each 
case in the matrix. Using analytical solutions, we 
were then able to evaluate the relative accuracy of 
each boundary condition. 

In the following section we review the mathe- 
matical formulation of the flat plate gust problem, 
and in Section III we summarize the numerical im- 
plementation. In Section IV we present and discuss 
the four FFBCs. Finally, in Section V we describe 
the details of our numerical study and present the 
results. 

II. Mathematical Formulation 

Governing Equation 

Consider an airfoil with chord length c in a flow 
with uniform upstream velocity in the ,r t di- 
rection. Let the fluid be an ideal gas which is in- 
viscid and non-heat-conducting. Far upstream, let 
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/I x = ac ,k {£ ~ ^ denote a small amplitude gust, 

where i is a unit vector in the X) direction. Here 
a = («.] where the amplitude |a| satisfies 

(a | << U,^, k = is the wave number vec- 

tor, and a and k satisfy a • k — 0 to ensure that the 
continuity equation is satisfied. 

We assume that the velocity field can be repre- 
sented by 


U(xJ) = U Q {x) + ii(xj) (2.1) 

where I o (x) is the steady potential velocity and 
u(x, t) is t he small, unsteady velocity. u(x, t) may 
be decomposed 9,10 into the sum of a known vortical 
component if R K and an unknown potential compo- 
nent V 0 , so that 


or • n. For a flat plate at zero degrees 

incidence, this reduces to 


00 

Ox 2 


—a <> e. 
or 




_ __ q ^ iff 1 (*£ I Ocv X ) + / A' 3 J.“;t 

Ox-2 


(2.G) 


for — £ < jq < §. and .r 2 — 0. 

Across the wake t he pressure is cont inuous, but 
0 has a discontinuity A 0 due to the unsteady cir- 
culation. Applying (2.4) on each side of (lie vortex 
sheet, behind the airfoil leads to 


Do 

Dt 


(A 0) = l) 


( 2 . 8 ) 


u(.r : t) = i? H) (2.2) 

The vortical component is essentially a func- 
tion of the mean flow Lagrangian coordinates and 
their spatial gradients. However, in the case of thin, 
uncambered airfoils at zero degrees incidence, the 
unsteady velocity decouples from the mean flow and 
ifl R 1 reduces to a x . 

The potential 0 satisfies the convective wave 
equation 

= (2.3) 

Dt c 0 2 Dt po po 

where is the convective derivative associated with 
the mean flow, and c*o and po are, respectively, the 
mean flow speed of sound and density. The unsteady 
pressure is given by 


P = 


~Po 0?) 


D o 0 
Dt ' 


(2.4) 


For non uni form mean flows, the coefficients of 
equation (2.3) and its source term depend strongly 
on the mean flow quantities. However, for a thin air- 
foil in a uniform parallel flow, the coefficients decou- 
ple from the mean flow and the source term vanishes. 
In this case the equation reduces to the constant, co- 
efficient convective wave equation 


1 A) 2 0 

co 2 nr 2 


- V 2 0 = 0. 


(2.5) 


Boundary Conditions 


for x\ > X 2 = 0. 

Finally, 0 must satisfy 

V0 — > 0 as .ri — ► — oc. (2.9) 


III. Numerical Scheme 


Let lengths be normalized by ^ , time by ju~^ 
0 by §|a|, a by |a|, and and c 0 by . Then 
equation (2.5) and boundary conditions (2.7) - (2.9) 
become 


d< P — _ n JfrlU’l -“f) + *A’3X3 

dx 2 — U2C 


for - 1 < r, < 1, x -2 = 0 


-^(A0) = 0 for xi > 1, x -2 = 0 


(3.1) 

(3.2) 
(3-3) 


Y(p — » 0 as .ri — > — oc (3.4) 

where = Jj 4- and M denotes the Mach num- 
ber. The normalized wave number k\ = ^ jf—, where 
v is the angular frequency of the gust, is called the 
reduced frequency. 

Introducing the change of variables 


;r = x\ , y — fit r 2 , £ = ;r 3 where fi 2 = 1 — M 2 

(3.5a) 

and 

- k-i M 2 

i’ = <)>ee*~ x e iklt - ik * z (3.5b) 

equation (3.1) is transformed into the frequency do- 
main and reduced to a Helmholtz equation 


At the airfoil surface the normal velocity com- 
ponent must vanish, so that (v^ Ri 4- Vo) ■ n = 0, 


d 2 tl' 
0 te* 


d 2 \b 

dy 2 


K 2 xb = 0 


(3.6a) 
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where 




Boundary conditions (3.2) - (3.4) become 


— - _ £2 x/d 2 for - 1 < .r < 1. // = 0, 

0y ti - - . 

(3.7) 

(”^rf + = 0 i »r .r > 1. // = U- 

p 2 ox 

(3.8) 

Yxp-tQ as .r — oo. (3.9) 

To solve (3.6) - (3.9) numerically, we introduce 


elliptic coordinates (//,£), whe e 


x = cos(7TJ7)cosli(7r^) 

(3.10a) 

11 

jn 

7 ? 

^3 

(3.10b) 

Equation (3.6) then becomes 


d 2 ii' 0 2 ij? i 

dif ‘ 1 

(3.11) 

7r 2 K 2 [sin* 2 (7r//) — sinh 2 (7C)]^ = 0 J 


and the airfoil boundary cond non (3.7) is 


— = -^ls\n(7Tn)e tkl (3.12) 
0 £ p 

for 0 < r) < L £ = 0. 

The wake condition (3.8) can be integrated to 
obtain 

At/? = Atbx e x > 1, (j = 0, 

(3.13) 

where the subscript t.e. denotes the airfoil trailing 
edge. For a flat plate, ?/; is an odd function of y. so 
that for points above and below the wake, ’& = ^ At/ ? 
and (3.13) becomes 

t p = ip t e e ikl{co * h{ "* ) ‘ l]/ t i * n = o, £ > 0. 

(3.14) 

Since ip is odd, it is necessary to solve (3.11) 
only in the upper half plane. One therefore imposes 
the symmetry condition ip = 0 on the streamline 
ahead of the airfoil. This coriesponds to the line q 
= 1,S>0. 

At the trailing edge, one must impose the Kutta 
condition. This requires that ^ = 0 at q — 0, £ = 
0. 

Finally, we note that condition (3.9) cannot be 
imposed at any finite distance from the airfoil with- 
out causing acoustic waves tn reflect off the outer 
grid boundary, thereby contaminating the numeri- 
cal solution. Instead, one must impose a radiation 


boundary condition which allows outgoing acoustic 
waves to exit the computational domain. We discuss 
this in the next section. 

Upon selection of an appropriate FFBO. the 
boundary value' problem can be discretized in a 
straightforward manner using finite differences. The 
resulting matrix equation can then be solved using 
a sparse matrix solver. See [4] for details. 


IV. Radiation Boundary Conditions 


The simplest radiation boundary condition is 
the Sommerfeld condition 11 . 



where u satisfies the Helmholtz equation 


d 2 u 

U7 2 


d' 2 u 

Oy 2 


k~ a 


0 


(4.1) 


(4.2) 


and the variable r denotes polar distance. Although 
this condition is exact at infinity, it leads to spurious 
reflections when imposed at a finite distance 12 . 

Many authors have developed improved radi- 
ation conditions to reduce the spurious reflections. 
Our purpose here is not to review the various condi- 
tions that have been developed, but to examine the 
performance of two advanced radiation conditions 
versus the Sommerfeld condition. 

We consider the Bayliss-Turkel' (BT) and 
Hagstrom-Hariharan 8 (HH) conditions of order one. 
These can be expressed conveniently in operator no- 
tation. as follows. First, let 


B = ik - — . 

2r 


Then the BT condition is 

0 


fe - s )" = 


0. 


Second, let 


H = ik 

r 


and 


“ ~{' k - u) 


H 


Then t he HH condition is 


( j - "I 


+ 


1 x 

2 r 2 Off 2 J 


8 r 2 


u = 0 


(4.3) 

(4.4) 
(4.5a) 

(4.5b) 

(4-6) 


where 0 denotes the polar angle. 

Conditions (4.1), (4.4), and (4.6) are accu- 
rate to 0(r _3/2 ), 0(r~ 5/2 ), and 0{r~ 7/2 ). respect- 
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ively'- 8 * 14 . It must be noted, however, that these 
conditions were derived for circular boundaries. The 
elliptic coordinate transformation (3.10) provides an 
outer boundary which is nearly circular for most 
problems. However, there will be a loss of accuracy 
as the outer boundary deviates from a circle. 

Each of the above conditions can be applied to 
the transformed potential since ij: satisfies the 
Helmholtz equation. However, it has been previ- 
ously observed that this does not lead to accurate 
far- field pressure calculations' 1,13 . A far superior ap- 
proach is to apply the radiation condition to the 
pressure p, which also satisfies equation (3.0). This 
leads to an additional operator acting on ?/; in the 
far-field boundary condition. From equation (2.4), 
this operator is 

(it- ~ A ) (4 - 7a) 


ik] 

4 = ~W' 


(4.7b) 


We may therefore identify the variable u in (4.1) - 
(4.0) with the pressure p and write (leaving out the 
minus sign) 


= &-*)*■ 


The Sommerfeld, BT, and HH pressure radiation 
conditions are then 






(4.10) 




respectively, where the frequency k in (4.1) - (4.0) is 
now identified with K as defined in (3.6b). 


V. Numerical Stuch 


The main objectives of our study were to deter- 
mine which radiation condition was most accurate, 
which was least sensitive to changes in grid, and 
which was best overall in terms of both accuracy and 
efficiency. We also wanted to compare the partial 
Sommerfeld (PS) condition, which uses only radial 
derivatives, with the full Sommerfeld (FS) condition, 
which retains all derivatives. 

To test the boundary conditions, we con- 
structed a test matrix of problems in which the Mach 


number M ranged from 0.1 to 0.5 to 0.8, with re- 
duced frequency values ki of 0.1, 0.5, 1.0, 2.0. and 
3.0 at each Mach number. For each case in the ma- 
trix we solved the gust problem on a set of 15 dif- 
ferent grids. The grids varied in their spacing in the 
tangential direction and also in the location of the 
outer grid boundary. The spacing in the radial di- 
rection was kept fixed at 24 points per wavelength. 
Each grid was designed to provide an accurate calcu- 
lation so that differences in accuracy that appeared 
would be due to the given FFBC. Two parameters in 
the GUST3D code were used to vary the grids. The 
parameter “neta”, which specifies the number of q 
grid points, was used to vary the spacing in the tan- 
gential direction. The parameter “mvaves”, which 
specifies the number of wavelengths t.o be used in 
the wake boundary condition, was used to vary the 
location of the outer grid boundary. Tables T and IT 
summarize the grid configurations that were used. 


Grid No. 

neta 

n waves 

1 

31 

3.0 

2 

31 

3.5 

3 

31 

4.0 

4 

31 

4.5 

5 

31 

5.0 

6 

36 

3.0 

i 

36 

3.5 

8 

36 

4.0 

9 

36 

4.5 

10 

36 

5.0 

11 

41 

3.0 

12 

41 

3.5 

13 

41 

4.0 

14 

41 

4.5 

15 

41 

5.0 


Table I Grid Parameters for ki =0.1, 0.5, 1.0 


Grid No. 

neta 

nwaves 

1 

41 

5.0 

2 

41 

5.5 

3 

41 

6.0 

4 

41 

6.5 

5 

41 

7.0 

6 

46 

5.0 

7 

46 

5.5 

8 

46 

6.0 

9 

46 

6.5 

10 

46 

7.0 

11 

51 

5.0 

12 

51 

5.5 

13 

51 

6.0 

14 

51 

6.5 

15 

51 

7.0 


Table II Grid Parameters for k\ = 2.0, 3.0 
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To assess accuracy, we compared numerical and 
exact pressure values on a semicircle of radius two 
chord lengths centered at the airfoil center, as shown 
in Figure 1. Figures 2 and 3 show comparisons be- 
tween numerical and exact results from typical cal- 
culations. 

To quantify the error, we used two error norms 
- absolute L y c and relative / ■. . The absolute L x 
norm was calculated by 


max „ 


p n — P n 
r comp ^ex * 


i = 0,1,.... 360 


where the semicircle was discretized by 301 uni- 
formly spaced points and P c \\ nX[) and are the com- 
puted and exact, pressure at point n. respectively. 
The relative L> error was calculated by 


. /v3GO (p n __ DriV , /v‘^0 ( p n )' 2 
Y ^n=0 \ r comp 'ex ) j V ' 

Using the above error norms, we were able to 
compare the accuracy of one FFBC versus another. 
We illustrate this in Figures 4-9. In Figure 4. for 
example, we show the absolute L oo error for all four 
boundary conditions for the case M = 0.K k\ = 1.0, 
on all 15 grids. From the ph i one can readily see 
which FFBC was most accurate and which one was 
least sensitive to the change in grid. At the top of 
the figure we show the maximum error, mean 
Lzc error, and standard deviation for each of the 
four boundary conditions. It is seen that for this 
particular case, the Hagstron Hariharan condition 
was superior to the other three conditions in all three 
categories. On the other hand, the Bayliss-Turkel 
condition had the smallest max and mean error in 
Figures 5 and 8, and the smallest standard deviation 
in Figure 9. 

Plots analogous to those in Figures 4-9 were 
made for every case in the test matrix, however for 
brevity they are omitted here This made it possi- 
ble to do a complete evaluation of the four bound- 
ary conditions. Tables III and IV summarize which 
FFBC was best in the categories of smallest max er- 
ror, smallest mean error, and smallest standard de- 
viation for all the cases in the t est matrix. It is clear 
that the Bayliss-Turkel and Hagstrom-Hariharan 
conditions performed the best. How much better 
they performed is shown in Figures 10 - 15. In Fig- 
ure 10, for example, we consider the max error. 
The figure shows the ratio of max error to the partial 
Sommerfeld max error for the other three boundary 
conditions. The ratio is shown for each case in the 
test matrix. An error ratio h ss than one indicates 
an improvement over the part al Sommerfeld condi- 
tion, which is the simplest o’ the four conditions. 
An inspection of Figures 10 - 15 indicates that the 


M 

k± 

Max Error 

Mean Error 

Stand Dev 

0.1 

0.1 

HH 

HH 

HH 

0.1 

0.5 

BT 

HH 

BT 

0.1 

1.0 

HH 

HH 

HH 

0.1 

2.0 

HH 

HH 

HH 

0.1 

3.0 

HH 

HH 

HH 

0.5 

0.1 

BT 

BT 

HH 

0.5 

0.5 

BT 

BT 

BT 

0.5 

1.0 

BT 

BT 

HH 

0.5 

2.0 

HH 

HH 

BT 

0.5 

3.0 

HH 

HH 

BT 

0.8 

0.1 

BT 

BT 

BT 

0.8 

0.5 

BT 

BT 

HH 

0.8 

1.0 

HH 

HH 

BT 

0.8 

2.0 

HH 

PS 

BT 

0.8 

3.0 

PS 

PS 

BT 

Table 

Based 

III Best Far-Field 

on Absolulte L x Error 

Boundary Condition 


Af 

k L 

Max Error 

Mean Error 

Stand Dev 

0.1 

0.1 

BT 

BT 

BT 

0.1 

0.5 

BT 

BT 

BT 

0.1 

1.0 

HH 

HH 

HH 

0.1 

2.0 

HH 

HH 

HH 

0.1 

3.0 

HH 

HH 

HH 

0.5 

0.1 

BT 

BT 

BT 

0.5 

0.5 

BT 

BT 

HH 

0.5 

1.0 

BT 

BT 

HH 

0.5 

2.0 

HH 

HH 

BT 

0.5 

3.0 

HH 

HH 

BT 

0.8 

0.1 

BT 

BT 

BT 

0.8 

0.5 

HH 

BT 

HH 

0.8 

1.0 

HH 

HH 

BT 

0.8 

2.0 

PS 

PS 

PS 

0.8 

3.0 

PS 

PS 

FS 


Table TV Best Far-Field Boundary Condition 
Based on Relative L 2 Error 

Bayliss-Turkel and Hagstrom-Hariharan conditions 
performed significantly better than the partial Som- 
rnerfeld condition both in terms of accuracy and in 
terms of consistency of solution (i.e., insensitivity 
to change in grid). On the other hand, the results 
also show that the part ial Sommerfeld condition per- 
formed as well as or better than the full Sommerfeld 
condition. 

To determine which FFBC was most accurate 
and which was least sensitive to changes in grid, we 
use the results shown in Tables III and IV. One finds 
that in the categories of lowest max error and lowest 
mean error, the Hagstrom-Hariharan condition was 
the best 29 out of GO times, while the Bayliss-Turkel 
condition was best 24 times. We therefore con- 
clude that the Hagstrom-Hariharan condition was 
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the most, accurate overall. On the other hand, in the 
category of smallest standard deviation, one finds 
that the Bayliss-Turkel condition was the best 15 
out of 30 times, while the Hagstrom-Hariharan con- 
dition was best 13 times. We t herefore conclude that 
the Bayliss-Turkel condition was the least sensitive 1 
overall to a change in grid. 

In terms of computational efficiency, partial 
Sornrnerfeld was the most efficient, of the four bound- 
ary conditions. This is due to its smaller discretiza- 
tion stencil, resulting in a sparser matrix than that 
of the other boundary conditions. The full Sornmer- 
feld and Bayliss-Turkel conditions required an in- 
crease in CPU time over part ial Sornrnerfeld of about 
3%. However, they required an increase in computer 
storage of about 52%. Similarly, the Hagstrom- 
Hariharan condition required an increase in CPU 
time over partial Sornrnerfeld of about 10%, while 
requiring an increase in computer storage of about 
102%. Clearly, the biggest cost associated with the 
three more advanced boundary conditions was the 
substantial increase in required storage. Since the 
Bayliss-Turkel condition was only slightly less accu- 
rate than the Hagstrom-Hariharan condition while 
providing the most grid independent solution at a 
cost of only a 52% increase in storage, we conclude 
that the Bayliss-Turkel condition was the best over- 
all boundary condition. 

Summary 

A detailed study of the accuracy and efficiency 
of the partial Sornrnerfeld, full Sornrnerfeld, order- 1 
Bayliss-Turkel, and order- 1 Hagstrom-Hariharan ra- 
diation boundary conditions was carried out. Each 
boundary condition was tested on a matrix of prob- 
lems with three different Mach numbers and five 
different reduced frequencies. It was determined 
that Hagstrom-Hariharan was most accurate over- 
all, Bayliss-Turkel was least, sensitive to changes in 
grid, and Bayliss-Turkel was best overall in terms 
of both accuracy and efficiency. In addition, it was 
determined that full Sornrnerfeld did not yield any 
improvements over partial Sornrnerfeld. Finally, it 
was determined that for far-field pressure calcula- 
tions, both Bayliss-Turkel and Hagstrom-Hariharan 
are sufficiently superior to Sornrnerfeld that they are 
well worth the additional cost of implementation. 
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Figu e 1 Semicircle of radius two chord lengths for 
far-field pressure calculation 



Figure 2 Comparison of far-field pressure on semicircle 
of radius two chord lengths. M = 0.5, kl = 1.0, 
flat plate in a transverse gust. 
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Figure 3 Comparison of far-field pressure on semicircle 
of radius two chord lengths. M = 0.5, kl = 3.0, 
flat plate in a transverse gust. 
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Figure 4 Absolute Ljnf error for all four FFBC’s on 1 5 
different grids for the case M=0.1, kl=1.0 
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Figure 5 Absolute L_inf error for all four FFBC's on 1 5 
different grids for the case M=0.5, kl = 1.0 
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Figure 6 Absolute L_inf error for all four FFBC’s on 15 
different grids for the case M=0.8, k 1=1.0 
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Relative L2 Error 
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Figure 7 Relative L_2 error for all four FFBCs on 15 
different grids for the case M=0. 1 , k 1 = 1 .0 
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Figure 8 Relative L_2 error for all four FFBC’s on 15 
different grids for the case M=0.5, kl=1.0 
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F gure 9 Relative L_2 error for all four FFBC's on 15 
different grids for the case M=0.8, k 1 = 1 .0 



Figure 10 Ratio of L_inf max error lo partial Sommerfeld 
LJnf max error for all cases in the test matrix 
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Figure 1 1 Ratio of L_inf mean error to partial Sommerfeld 
L_inf mean error for all cases in the test matrix 



Figure 12 Ratio of L_inf stand, dev. to partial Sommerfeld 
L_inf stand, dev. for all cases in the test matrix 
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Figure 13 Ratio of L_2 max error to partial Sommerfeld 
L_2 max error for all cases in the test matrix 



Figure 14 Ratio of L_2 mean error to partial Sommerfeld 
L_2 mean error for all cases in the test matrix 
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Figure 15 Ratio of L_2 stand, dev. to partial Sommerfeld 
L_2 stand, dev. for all cases in the test matrix 
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